
#include "MaxPoolingRBM.h"

// [*] nw_rbm::CUnitPooling

nw_rbm::CUnitPooling::CUnitPooling()
{
	Clear();
}

void nw_rbm::CUnitPooling::Clear()
{
	State = -1;
	BiasValue = 0.0f;
	ExtraEnergy = 0.0f;
	Expectation = 0.5;
	//DerivativeFactor_SpaCod = 1.0;
}

// [*] nw_rbm::CRBM_with_MaxPooling

nw_rbm::CRBM_with_MaxPooling::CRBM_with_MaxPooling()
{
	ClearUnitSettings();
}

void nw_rbm::CRBM_with_MaxPooling::ClearUnitSettings()
{
	CRestrictedBoltzmannMachine::ClearUnitSettings();
	AmountUnit_Pooling = 0;
	ListUnit_Pooling.Clear();
	ListIndex_MP_Unit_PerHidden.Clear();
	ListIndex_MP_Brother_PerHidden.Clear();
}

void nw_rbm::CRBM_with_MaxPooling::SetState_of_Units_Pooling( const int* SrcPtrState )
{
	int i;
	CUnitPooling* pU = ListUnit_Pooling.pElement;
	for ( i = 0; i < AmountUnit_Pooling; ++i )
		pU[i].State = SrcPtrState[i];
}

void nw_rbm::CRBM_with_MaxPooling::SetState_of_Units_Pooling( const CSimpleTypeArray<int>& SrcListState )
{
	int i;
	CUnitPooling* pU = ListUnit_Pooling.pElement;
	for ( i = 0; i < AmountUnit_Visible; ++i )
		pU[i].State = SrcListState.pElement[i];
}

void nw_rbm::CRBM_with_MaxPooling::GetState_of_Units_Pooling( RBM_FLOAT* DstPtrState ) const
{
	int i;
	const CUnitPooling* qU = ListUnit_Pooling.pElement;
	for ( i = 0; i < AmountUnit_Pooling; ++i )
		DstPtrState[i] = (qU[i].State>=0 ? 1.0f : 0.0f);
}

const RBM_FLOAT* nw_rbm::CRBM_with_MaxPooling::GetState_of_Units_Pooling( CSimpleTypeArray<RBM_FLOAT>& DstListState ) const
{
	DstListState.Resize(AmountUnit_Pooling);
	int i;
	const CUnitPooling* qU = ListUnit_Pooling.pElement;
	for ( i = 0; i < AmountUnit_Pooling; ++i )
		DstListState.pElement[i] = (qU[i].State>=0 ? 1.0f : 0.0f);
	return DstListState.pElement;
}

void nw_rbm::CRBM_with_MaxPooling::CheckPoolingState_and_UpdateHiddenState()
{
	int i, ee, j;
	const int* qUI;
	const int* qWI;
	int nU;
	CUnitBasic* pUnit_Hidden = ListUnit_Hidden.pElement;

	for ( i = 0; i < AmountUnit_Pooling; ++i )
	{
		CUnitPooling& tUnit_Pooling = ListUnit_Pooling.pElement[i];
		nU = tUnit_Pooling.Dimension;
		qUI = tUnit_Pooling.vTargetHandle.pElement;
		qWI = tUnit_Pooling.vEdgeHandle.pElement;
		if (tUnit_Pooling.State >= nU) tUnit_Pooling.State = -1;
		for ( ee = 0; ee < nU; ++ee )
		{
			j = qUI[ee];
			pUnit_Hidden[j].State = (ee==tUnit_Pooling.State ? 1.0f : 0.0f);
		}
	}
}

void nw_rbm::CRBM_with_MaxPooling::InitAmountUnit_Pooling( int AmountUnit_Pooling )
{
	this->AmountUnit_Pooling = AmountUnit_Pooling;
	CUnitPooling tUnit;
	tUnit.Clear();
	ListUnit_Pooling.Clear();
	ListUnit_Pooling.Resize(AmountUnit_Pooling, tUnit);
	ListIndex_MP_Unit_PerHidden.Clear();
	ListIndex_MP_Unit_PerHidden.Resize(AmountUnit_Hidden, -1);
	ListIndex_MP_Brother_PerHidden.Clear();
	ListIndex_MP_Brother_PerHidden.Resize(AmountUnit_Hidden, -1);
}

int nw_rbm::CRBM_with_MaxPooling::NewUnit_Pooling()
{
	int ret = AmountUnit_Pooling++;
	CUnitPooling tUnit;
	tUnit.State = -1;
	ListUnit_Pooling.PushBack(tUnit);
	return ret;
}

int nw_rbm::CRBM_with_MaxPooling::NewUnit_Hidden(
	bool fIsRealValued, RBM_FLOAT Sigma, int WeightItemIndex_Hidden )
{
	fIsRealValued = false;
	Sigma = 1.0f;
	ListIndex_MP_Unit_PerHidden.PushBack(-1);
	ListIndex_MP_Brother_PerHidden.PushBack(-1);
	return CRestrictedBoltzmannMachine::NewUnit_Hidden(fIsRealValued, Sigma, WeightItemIndex_Hidden);
}

int nw_rbm::CRBM_with_MaxPooling::NewUnit_Hidden(
	int WeightItemIndex_Hidden )
{
	ListIndex_MP_Unit_PerHidden.PushBack(-1);
	ListIndex_MP_Brother_PerHidden.PushBack(-1);
	return CRestrictedBoltzmannMachine::NewUnit_Hidden(false, 1.0f, WeightItemIndex_Hidden);
}

void nw_rbm::CRBM_with_MaxPooling::NewRelation_MaxPooling( int UnitIndex_Hidden, int UnitIndex_Pooling )
{
	if (UnitIndex_Hidden<0 || UnitIndex_Hidden>=AmountUnit_Hidden) return;
	if (UnitIndex_Pooling<0 || UnitIndex_Pooling>=AmountUnit_Pooling) return;
	CUnitPooling& tUnitPooling = ListUnit_Pooling.pElement[UnitIndex_Pooling];
	int* pIMP_U = ListIndex_MP_Unit_PerHidden.pElement + UnitIndex_Hidden;
	int* pIMP_B = ListIndex_MP_Brother_PerHidden.pElement + UnitIndex_Hidden;
	if (*pIMP_U >= 0) return;
	*pIMP_U = UnitIndex_Pooling;
	*pIMP_B = tUnitPooling.Dimension;
	tUnitPooling.NewConnection(UnitIndex_Hidden, -1);
}

const RBM_FLOAT* nw_rbm::CRBM_with_MaxPooling::MakeHiddenDerivative_for_SparseCoding(
	CSimpleTypeArray<RBM_FLOAT>& vHiddenDerivative )
{
	int j;
	vHiddenDerivative.Resize(AmountUnit_Hidden);
	RBM_FLOAT* ret = vHiddenDerivative.pElement;
	const int* qPI = ListIndex_MP_Unit_PerHidden.pElement;
	const CUnitPooling* qUP = ListUnit_Pooling.pElement;
	for ( j = 0; j < AmountUnit_Hidden; ++j )
	{
		CUnitBasic& tHj = ListUnit_Hidden.pElement[j];
		if (qPI[j] < 0)
			ret[j] = 2.0f * tHj.Expectation * (1.0f - tHj.Expectation);
		else
			ret[j] = 2.0f * tHj.Expectation * (1.0f - qUP[qPI[j]].Expectation);
			//ret[j] = 2.0f * tHj.Expectation * qUP[qPI[j]].DerivativeFactor_SpaCod;
	}
	return ret;
}

void nw_rbm::CRBM_with_MaxPooling::UpdateUnitsBiasValue_Pooling()
{
	int i, ee, j;
	double ss;
	const int* qUI;
	const int* qWI;
	int nU;
	CUnitBasic* pUnit_Hidden = ListUnit_Hidden.pElement;

	for ( i = 0; i < AmountUnit_Pooling; ++i )
	{
		CUnitPooling& tUnit_Pooling = ListUnit_Pooling.pElement[i];
		nU = tUnit_Pooling.Dimension;
		qUI = tUnit_Pooling.vTargetHandle.pElement;
		qWI = tUnit_Pooling.vEdgeHandle.pElement;
		ss = -MAXIMUM_DOUBLE_NUM_V;
		for ( ee = 0; ee < nU; ++ee )
		{
			j = qUI[ee];
			ss = CMyMath::log_Sum(ss,
				this->GetUnitBiasValue_Hidden(j));
		}
		tUnit_Pooling.BiasValue = (RBM_FLOAT)ss;
	}
}

void nw_rbm::CRBM_with_MaxPooling::ClearMaxPoolingUnitsExtraEnergy()
{
	int i;
	for ( i = 0; i < AmountUnit_Pooling; ++i )
	{
		CUnitPooling& tUnit_Pooling = ListUnit_Pooling.pElement[i];
		tUnit_Pooling.ExtraEnergy = 0.0f;
	}
}

void nw_rbm::CRBM_with_MaxPooling::ComputeExpectation_of_Units_Hidden()
{
	int i, ee, j;
	double ss, tt, pp;
	const int* qUI;
	const int* qWI;
	int nU;
	CUnitBasic* pUnit_Hidden = ListUnit_Hidden.pElement;
	const int* qIMP_U = ListIndex_MP_Unit_PerHidden.pElement;
	const int* qIMP_B = ListIndex_MP_Brother_PerHidden.pElement;

	for ( j = 0; j < AmountUnit_Hidden; ++j )
	{
		if (qIMP_U[j] < 0)
			pUnit_Hidden[j].ComputeExpectation();
	}

	for ( i = 0; i < AmountUnit_Pooling; ++i )
	{
		CUnitPooling& tUnit_Pooling = ListUnit_Pooling.pElement[i];
		nU = tUnit_Pooling.Dimension;
		qUI = tUnit_Pooling.vTargetHandle.pElement;
		qWI = tUnit_Pooling.vEdgeHandle.pElement;
		ss = 0.0f - tUnit_Pooling.ExtraEnergy;
		for ( ee = 0; ee < nU; ++ee )
		{
			j = qUI[ee];
			ss = CMyMath::log_Sum(ss, pUnit_Hidden[j].Energy);
		}
		//tUnit_Pooling.DerivativeFactor_SpaCod = (RBM_FLOAT)CMyMath::my_exp(0.0f - ss);
		pp = 0.0f;
		for ( ee = 0; ee < nU; ++ee )
		{
			j = qUI[ee];
			tt = CMyMath::my_exp(CMyMath::log_Product(pUnit_Hidden[j].Energy, 0.0f-ss));
			pUnit_Hidden[j].Expectation = (RBM_FLOAT)tt;
			pp += tt;
		}
		tUnit_Pooling.Expectation = (RBM_FLOAT)pp;
	}
}

void nw_rbm::CRBM_with_MaxPooling::ComputeBiasOnlyExpectation_of_Units_Hidden()
{
	CSimpleTypeArray<RBM_FLOAT> vHiddenEnergy;
	const RBM_FLOAT* qHE = this->GetEnergy_of_Units_Hidden(vHiddenEnergy);
	this->SetBiasOnlyEnergy_of_Units_Hidden();
	this->ComputeExpectation_of_Units_Hidden();
	this->SetEnergy_of_Units_Hidden(vHiddenEnergy);
}

void nw_rbm::CRBM_with_MaxPooling::SampleState_of_Units_Hidden()
{
	int i; // , j;
	RBM_FLOAT pp;
	CUnitBasic* pUnit_Hidden = ListUnit_Hidden.pElement;
	const int* qIMP_U = ListIndex_MP_Unit_PerHidden.pElement;
	const int* qIMP_B = ListIndex_MP_Brother_PerHidden.pElement;

	//for ( j = 0; j < AmountUnit_Hidden; ++j )
	//{
	//	if (qIMP_U[j] < 0)
	//		pUnit_Hidden[j].SampleState();
	//}

	for ( i = 0; i < AmountUnit_Pooling; ++i )
	{
		CUnitPooling& tUnit_Pooling = ListUnit_Pooling.pElement[i];
		pp = (RBM_FLOAT)CMyRandom::NextDouble();
		tUnit_Pooling.State = (pp<tUnit_Pooling.Expectation ? 0 : -1);
	}

	SampleState_of_Units_Hidden_GivenPoolingState();
}

void nw_rbm::CRBM_with_MaxPooling::SampleState_of_Units_Hidden_GivenPoolingState()
{
	int i, ee, j;
	RBM_FLOAT pp, tt;
	const int* qUI;
	const int* qWI;
	int nU;
	CUnitBasic* pUnit_Hidden = ListUnit_Hidden.pElement;
	const int* qIMP_U = ListIndex_MP_Unit_PerHidden.pElement;
	const int* qIMP_B = ListIndex_MP_Brother_PerHidden.pElement;

	for ( j = 0; j < AmountUnit_Hidden; ++j )
	{
		if (qIMP_U[j] < 0)
			pUnit_Hidden[j].SampleState();
	}

	for ( i = 0; i < AmountUnit_Pooling; ++i )
	{
		CUnitPooling& tUnit_Pooling = ListUnit_Pooling.pElement[i];
		nU = tUnit_Pooling.Dimension;
		qUI = tUnit_Pooling.vTargetHandle.pElement;
		qWI = tUnit_Pooling.vEdgeHandle.pElement;
		if (tUnit_Pooling.State < 0)
		{
			tUnit_Pooling.State = -1;
			for ( ee = 0; ee < nU; ++ee )
			{
				j = qUI[ee];
				pUnit_Hidden[j].State = 0.0f;
			}
		}
		else
		{
			pp = (RBM_FLOAT)CMyRandom::NextDouble() * tUnit_Pooling.Expectation;
			tUnit_Pooling.State = -1;
			for ( ee = 0; ee < nU; ++ee )
			{
				j = qUI[ee];
				tt = pUnit_Hidden[j].Expectation;
				if (tUnit_Pooling.State<0 && (pp<tt || ee==nU-1))
				{
					tUnit_Pooling.State = j;
					pUnit_Hidden[j].State = 1.0f;
				}
				else
				{
					pp -= tt;
					pUnit_Hidden[j].State = 0.0f;
				}
			}
		}
	}
}

using namespace nw_rbm;

void TestMP_RBM()
{
	//CMyRandom::InitRandomSeed();

	char strTemp[10001] = "";

	ofstream outFile_GlobalLog("Global.log");
	CExternalTools::ConsoleTools::MakeCurrentLine("", true);
	CExternalTools::ConsoleTools::OutputTextToStream(outFile_GlobalLog, "", true);
	ofstream outFile_ShortLog("Short.log");

	int NV = 4;
	int NH = 4;
	int AmountIteration = 550;
	RBM_FLOAT LearningRate = 0.3f;
	RBM_FLOAT DecayRate = 0.0f;
	RBM_FLOAT Momentum = 0.5f;
	RBM_FLOAT SparseLevel = 0.02f;
	RBM_FLOAT SparseCodingLambda = 3.0f;
	int AmountExample = 10000;

	//CRestrictedBoltzmannMachine RBM;
	CRBM_with_MaxPooling RBM;

	int i, j, ow, ov, oh;
	int i_iteration, i_example;
	RBM_FLOAT tet;
	RBM_FLOAT tVisErr, tHidPen;
	RBM_FLOAT VisibleError, HiddenPenalty;
	RBM_FLOAT DeltaConDiv, DeltaSpaCod;

	RBM.Clear_WeiUniTra_Settings();

	RBM.InitAmountUnit_Pooling(1);

	for ( i = 0; i < NV; ++i )
	{
		ov = RBM.NewWeightItem_Visible();
		RBM.NewUnit_Visible(true, 0.2f, ov);
	}
	for ( j = 0; j < NH; ++j )
	{
		oh = RBM.NewWeightItem_Hidden();
		RBM.NewUnit_Hidden(false, 1.0f, oh);
		RBM.NewRelation_MaxPooling(j, 0);
	}

	RBM.SetConnectionManagers(
		new CConnectionManager_Basic,
		new CConnectionManager_Basic);

	for ( j = 0; j < NH; ++j ) for ( i = 0; i < NV; ++i )
	{
		ow = RBM.NewWeightItem_Edge();
		RBM.NewRelation_Edge(i, j, ow);
	}
	RBM.PrintEdgeWeightAmountPerLine = NV;

	RBM.InitWeightItems(0.1f);

	CSimpleTypeArray<RBM_FLOAT> ExampleTypeList;
	ExampleTypeList.Resize(AmountExample);
	for ( i_example = 0; i_example < AmountExample; ++i_example )
		ExampleTypeList.pElement[i_example] = (RBM_FLOAT)CMyRandom::NextDouble();

	CSimpleTypeArray<RBM_FLOAT> InputList;
	InputList.Resize(NV);

	DeltaConDiv = DeltaSpaCod = 0.0f;

	for ( i_iteration = 0; i_iteration < AmountIteration; ++i_iteration )
	{
		VisibleError = 0.0f;
		HiddenPenalty = 0.0f;
		for ( i_example = 0; i_example < AmountExample; ++i_example )
		{
			InputList.Fill(0.0f);
			tet = ExampleTypeList.pElement[i_example];
			if (tet < 0.333)
			{
				InputList.pElement[0] = 1.0f;
				InputList.pElement[2] = 1.0f;
			}
			else if (tet < 0.666)
			{
				InputList.pElement[1] = 1.0f;
				InputList.pElement[3] = 1.0f;
			}

			RBM.SetState_of_Units_Visible(InputList);
			tVisErr = RBM.ComputeGradient_ConDiv_GivenUnitsVisible(NULL, false, 1.0f, 1.0f, 1);
			RBM.SetState_of_Units_Visible(InputList);
			tHidPen = RBM.ComputeGradient_SpaCod_GivenUnitsVisible(NULL, SparseLevel, 1.0f, 1.0f);

		}
		VisibleError /= AmountExample;
		HiddenPenalty /= AmountExample;

		DeltaConDiv = RBM.Update_ContrastiveDivergence(
			LearningRate, DecayRate, DecayRate, Momentum);
		DeltaSpaCod = RBM.Update_SparseCoding(
			LearningRate, SparseCodingLambda);

		sprintf_s(strTemp, "%d\t%.4lf\t%.4lf\t%.4lf\t%.4lf\t%.4lf", i_iteration+1,
			DeltaConDiv, VisibleError, DeltaSpaCod, HiddenPenalty,
			VisibleError + HiddenPenalty * SparseCodingLambda);

		if (((i_iteration+1)^i_iteration)==(i_iteration+i_iteration+1))
		{
			CExternalTools::ConsoleTools::OutputTextToStream(outFile_GlobalLog, strTemp, true);
			DeltaConDiv = DeltaConDiv;
		}
		if (i_iteration==AmountIteration-1 || CMyRandom::NextInt(AmountIteration/200)==0)
		{
			CExternalTools::ConsoleTools::MakeCurrentLine(strTemp, false);
			DeltaConDiv = DeltaConDiv;
		}
	}
	RBM.SaveWeightModelToStream(outFile_ShortLog);

	CExternalTools::ConsoleTools::OutputTextToStream(outFile_GlobalLog, "", true);

	sprintf_s(strTemp, "%.6lf", VisibleError + HiddenPenalty * SparseCodingLambda);
	CExternalTools::ConsoleTools::OutputTextToStream(outFile_GlobalLog, strTemp, true);
	CExternalTools::ConsoleTools::MakeCurrentLine(strTemp, true);

	CExternalTools::ConsoleTools::MakeCurrentLine("", true);

	outFile_ShortLog.close();
	outFile_ShortLog.clear();
	CExternalTools::ConsoleTools::MakeCurrentLine("", true);
	CExternalTools::ConsoleTools::OutputTextToStream(outFile_GlobalLog, "", true);
	outFile_GlobalLog.close();
	outFile_GlobalLog.clear();
}

